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Abstract. Cosmic microwave background observations are most commonly 
analyzed by estimating the power spectrum. In the limit where the CMB statis- 
tics are perfectly Gaussian, this extracts all the information, but the CMB also 
contains detectable non-Gaussian contributions from secondary, and possibly 
primordial, sources. We review possible sources of CMB non-Gaussianity and 
describe statistical techniques which are optimized for measuring them, comple- 
menting the power spectrum analysis. The machinery of A'^-point correlation 
functions provides a unifying framework for optimal estimation of primordial 
non-Gaussian signals or gravitational lensing. We review recent results from 
applying these estimators to data from the WMAP satellite mission. 



1. Introduction 

Observations of the cosmic microwave background have had an enormous im- 
pact on our understanding of cosmology and the level of precision with which 
cosmological parameters can be constrained. B eginning with first detection of 
anisotropy on large angular scales by COBE ( Smoot et al.l 1199^ ). successive 



generations of experiments have provided increasingly p recise measurements 



of CMB temperature f iucutations on small e r scales (e . g . Hananv et al.| 1 20001 : 
Netterfield et al. 2002; H inshaw et aD l2003l : iKuo et al.1 l2007l : lAde et al.l boosl : 
Lueker et al. 2009. ; .Fowler et al. 2010l ). In particular, the WMAP satellite has 



measured CMB fluctuations on the full sky with high signal-to-noise for all an- 
gular scales larger than ~20 arcmin. 

It is convenient to represent the CMB in harmonic space. If r(n) denotes 
the value of the temperature in line-of-sight direction n, then we expand in 
spherical harmonics: 

oo e 

^(n) = E E aimYimin) (1) 

i=2m=-E 

to define the multipoles aim- In this representation, the most general rotationally 
invariant two-point correlation function is: 

(«£imi«£2m2) = { — ^)"^^ C£-^d£^£2^rni,-m2 (2) 

where this equation defines the power spectrum C^. 

CMB observations are typically analyzed via the power spectrum: the main 
result of the analysis is a measurement of Ci and an estimate of the statistical 
uncertainty, which can be compared with theoretical predictions for the power 
spectrum. This approach is powerful because the CMB is a Gaussian field to a 
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good approximation. Gaussianity implies that the full probability distribution 
for the CMB map (i.e. the multivariate PDF for the variables a^m) is determined 
by the two-point correlation function, so that the power spectrum contains all 
the information. However, there are weak sources of non-Gaussianity which 
contain qualitatively new cosmological information, and require statistical tech- 
niques which complement the power spectrum. In this article, we will review 
sources of non-Gaussianity either in the early universe during inflation, or in 
the late universe via gravitational lensing. Our emphasis will be on statistical 
techniques for extracting these non-Gaussian signals, and we will present results 
from applying these techniques to WMAP data when possible. 

2. Primordial non-Gaussianity 

2.1. Three-point signals from inflation 

Consider the simplest model of inflation: single-field slow-roll inflation with 
standard kinetic term, 

c = ^R-l{d^<p){d,<p)-v{<^) (3) 

In such a model, the adiabatic curvature perturbation C(k) generated during 
inflation is a Gaussian field; in particular the three-point correlation function is 
zeroQ However, there are infiationary models which can generate a three-point 
function which is observationally distinguishable from zero. In this section we 
will review the phenomenology of such models. 

Let us first make some general mathematical comments. If we restrict 
attention to models which satisfy statistical translation and rotation invariance, 
then the three-point function of the curvature perturbation ( must be of the 
form 

(C(ki)C(k2)C(k3)) = F{h,k2, h) {27rf6^{k, + k2 + ks) (4) 

where the function F is called the bispectrum, and only depends on the lengths ki 
of the Fourier wavenumbers kj which form a closed triangle. For the infiationary 
models we will consider in this section, the bispectrum will have a simple de- 
pendence on overall scale (roughly F{aki,ak2,aks) ~ a~^~^'^^^~^''^ F{ki,k2,k3)) 
so that we can think of the bispectrum as a function of a 2-parameter family of 
triangle shapes without regard to overall scale. We will use the term "squeezed" 
for a triangle shape with ki <^ ^2 (which implies ki <C k^), "flattened" for a 
shape with ki ~ (/c2 + ^3)1 and "equilateral" for a shape with ki k. k2 ~ k^. 

The curvaton model is a multifleld model of inflation in which the source 
of primordial curvature fiuctuations is not the inflaton, but a second field a 
which does not dominate the energy density during inflation, but decays after 
the inflaton, when the cur v aton is oscillating near th e minimum of its potential 
(|Linde fc Mukhanovl ITqQTI : iMoroi k Takahashil [2OO1I : iLvth. Ungarelh fc WancS 



^This statement is an approximation, but a good one: the three-point functio n in such models 
corresponds roughly to ~ (5/12)(l — n^), far too small to be detectable (lAcguaviva et all 
[200g ; l^faIdacena^[200a ). 
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120031 1. The power spectrum Pf^Sa/a) of the fractional curvaton perturbation 
[5a /a) and the fraction / = p^/ Ptot of the energy density due to the curva- 
ton (both quantities evaluated at curvaton decay) are free parameters of the 
model. The perturbation to the curvaton energy density is given by [bpaj Pa) = 
2{5a/a) + {Sa/a)"^ , and can therefore be a non-Gaussian field if the quadratic 
term is non-negligible compared to the linear term. There is a region of pa- 
rameter space (i.e. / ^ 1 with Ps^/a = 0(10^^/^^)) in which the primordial 
curvature fluctuation (after curvaton decay) is significantly non-Gaussian and 
has a power spectrum which is consistent with CMB observations. The non- 
Gaussian curvature fiuctuation in the curvaton model can be written: 

C(x) = Cg(x) - \fSf\G{^f (5) 

where C,g is a Gaussian field, and /j^lf"' is a free parameter. (The factor 3/5 
is purely conventional.) If we write this equation in Fourier spac^, we get the 
following form of the bispectrum F: 

Fih,k2, h) = f§tF\oc.i{kiM. k3) (6) 

where 

6A? / 1 1 1 \ 

Flocdh,k2,k3) = ^ I pp + pp + pp (7) 

The local bispectrum is largest for squeez ed triangles. Conversely, there is a theo- 
rem, the single-field consistency relation ( Maldacenal 20031 : ICreminelli &: Zaldarriaga 



20041 ). which states that the bispectrum is always small (i.e. 0(10 ^)) in squeezed 



triangles, provided that inflation is single-field (but allowing for aribtrary self- 
interactions of the inflaton) . Thus the presence of a second field in the curvaton 
model is actually a necessary ingredient for any model which generates /Jvl^^ 
larger than 0(10"^). 

For single-field infiation, there is also a theorem ( Senatore. Smith &: Zaldarriagal 



2010l ) which shows that the most general primordial bispectrum is: 

F{ki,k2, ks) = /^Y'^equil(^l, ^2, ^s) + /^L^°^^orthog (^1 , ^2 , ^3) (8) 

where the "equilateral" and "orthogonal" bispectra are defined by|f| 

18A2 ^kl +k2- k3)ikl -k2 + k3){-kl + k2 + ks) 



%nil{kl,k2,k3) 



//^3 z^3 /-.3 



36A^ 1 



Forthogikl,k2,k3) = 3FEqui\ikl,k2,k3) -\ — ^ 1.2^21,2 (9) 





^Not e that in Eqs. (O and ([O]). we have assumed scale invariance, so that P(;{k) — A^k ^. 
See ISenatore. Smith fc Zaldarriagal l|2010f ) for expressions for the primordial bispectra with a 
power law spectrum P(;{k) = Ajfc~''(fe/feo)"''~^- 

''The equilateral and orthogonal bispectra defined in Eq. Q are actually approximations to 
more precise expressions for the bispectra generated during inflation. The approximation is 
made so that F{ki,k2,kz) will be a sum of a small number of terms which are factorizable 
in the form f{k\)g{k2)h{ki). As will be discussed in the next subsection, this factorizability 
condition is necessary in order to make the data analysis computationally tractable. 
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The equilateral shape is largest for equilateral triangles (as the name suggests) 
while the orthogonal shape changes sign between equilateral and flattened tri- 
angles. Both single-field shapes vanish in the squeezed limit, as required by the 
consistency relation. 



For a particular single-field model, such as DBI inflation (jAlishahiha. Silverstein Sz Tone 



I2OO4I ) or ghost inflation ( Arkani-Hamed et al. 20041 ) . the coefficients /^^'' and 



fNL^°^ in Eq. ([8]) can be calculated in terms of fundamental parameters of the 
model. In this article, our emphasis will be on data analysis. We simply remark 
that the forms of the bispectra Fiocai) -^equih -Porthog given above provide a point 
of contact between theory and data. Given data from a CMB experiment like 
WMAP, our job is to determine observational limits on the three fj\fL param- 
eters with statistical errors which are as small as possible. Once this has been 
done, any inflationary model can be compared with observations by calculating 
the fisiL parameters. 

Primordial non-Gaussianity is a particularly interesting probe of inflation 
because it can rule out qualitative classes of models. A robust detection of 
/i^cai _^ would rule out all single- fleld models. Detection of any nonzero Jnl 
parameter would rule out the simplest inflationary model, single-field inflation 
with standard kinetic term and slow-roll potential V{(f)). In the next two sub- 
sections we will turn to data analysis: we will construct optimal estimators for 
constraining the /nl parameters from data, and present results from WMAP. 

2.2. Estimators for primordial non-Gaussianity 

Let us begin by making some general statements about how three-point signals 
are estimated from data, before specializing to the case of the three inflationary 
shapes from the preceding subsection. 

We will make the approximation that the evolution from the 3D intitial 
curvature fluctuation C(k) to the 2D CMB aim is ^ linear operation, so that the 
three-point function of C translates linearly into the three-point function of the 
CMB: 

(10) 

where B^^^ 

,km2,i3m3 different for each X S {local, equil, orthog}. 
When we measure the CMB in a real experiment, we measure the a^m's 
plus some instrumental noise: 

obs true i noise /"I 1 ^ 

We will assume that the instrumental noise is Gaussian, and that the total 
covariance matrix of the observed a^^'s is given by the sum of a signal term 
(which is diagonal in i,m) and a noise term (which is non-diagonal): 

{O'timiO. 121x12) — (~l)™'^C'^i<J^i^2'^nii,-m2 ^i-i_mi/2m2 

(12) 

The detailed form of the noise covariance matrix Ni^rai,t2m2 encapsulates details 
of the experiment such as spatial gradients in instrumental noise level (repre- 
sented by assigning different noise variance to different pixels), the sky mask 
imposed to remove regions of high foreground emission (represented by assign- 
ing infinite noise variance to pixels which are masked), and the shapes of the 
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beams. We will denote the total covariance matrix on the RHS of Eq. (|12p by 
Ceimi,e2m2 ^iid its inverse by (C'-i)^i»"i/2m2^ fjij^^ matrix C represents all the 
properties of the experiment that we will need to know about in order to write 
down an estimator for /nl- 

Now we can ask the following general question: Given Bf^^^ ^^^^ ^^^^ and 

C^imi,^2m2) how should we estimate the value of the coefficient /j^^. ™ Eq. ([TO]) 
from the noisy CMB atml The optimal (i.e. mininum variance) estimator was 
found in ( Creminelli et al. 20061 ): 



f^L = J {t[C-M - {C-'ar^ (den.T[C~'a'])J (13) 
where the quantities T and dT are defined by ( Smith &: Zaldarriaga 20061 ): 



'^iC-'a] = ^J2BM2rn2,i,rnAC-'aY^"'^{C-'aY'"''iC-'aY^"'-^ 

li'mi 

d,mT[C-\] = ]^Y.BlmAm,/2m2{C-^^Y^"''{C-^<^Y^'"' (14) 

The normalization constant F, which is included so that f^i^ will be an unbiased 
estimator of /jy^ (i.e. {f^i) = /jvl) '^^^ determined by Monte Carlo. The 
estimator f^^ in Eq. (jl3p is a sum of a three-point term and a one-point "coun- 
terterm" which reduces the estimator variance in the presence of inhomogeneous 
noise. The one-point term vanishes if the noise is homogeneous, since rotation 
invariance implies {d£mT[C~^a']) = for i > 0. We will see more examples 
of this estimator structure, with an A^-point leading term plus counterterms of 
lower order, in ^3.3.1 

To evaluate the estimator f^i, we need algorithms for computing aim — > 
(C^^a)^™ and a — )■ T[a\. The first algorithm depends on details of the exper- 
iment being analyzed (via C~^) but not the form of the three-point function 
Be^mi,e2m2,i3m3'i we will discuss its implementation in WMAP in the next sub- 
section. The rest of this subsection is devoted to briefly describing the second 
algorithm, which depends only on the three-point function. 

For a completely general three-point function -B^imi,^2m2,^3m3) there is no 
algorithm for the operation a — )• T[a\ which is faster than the ©(^max) harmonic- 
space sum in Eq. (|14p . For current CMB experiments (with ima.x = O(IO^)) this 
is computationally prohibitive. However, if the three-point function is a sum of 
factorizable terms, in the sense that 

N 

R — z^^*) (^ 

-D£imi,£2m2,£3m3 — "^i^2^3^^1^2^3")'^1^2^3 



X < 



'(2£i + l) (2^2 + 1) (2^3 + 1) / h h h 



An V "^1 "^2 rns 



and N is not too la rge, then there is a fast (i.e. 0(£l,„ N)) algorithm for the 
a T[a\ operation ( Komatsu. Spergel Wandeltl l2005f : ICreminelh et al.ll2006l : 
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Smith &: Zaldarriaga 20061 ) . (A closely related algorithm exists for the operation 



a — >■ dirnT[a], which is also needed in order to evaluate the one-point term in 
the optimal estimator in Eq. (jl3p .) In practice, this factorizability condition 
is satisfied for "interesting" forms of the three-point function such as the local, 
equilateral and orthogonal shapes. 

In a little more detail, we can also define a 3D factorizability condition for 
the primordial bispectrum F{ki,k2,k3): we say that F is factorizable if it is a 
sum of terms 

N 

F{ki,k2, ks) = a,iki)i3,ik2hi{k3) (16) 
1=1 

It can be shown that a primoridal bispectrum F which is factorizable (in the 
sense defined by Eq. ()16p ) evolves to a CMB three-point function which is also 
factorizable (in the sense defined by Eq. (|15p ). although the number of terms 
is typically increased by a factor lO-lOOo 

The optimal estimator can also be used to search for three-point signals 
which are generated by secon dary sources of anisotropy s uch as point sources 
or gra v itational lensing (e.g. 'Spergel &: Goldbere? (l9W); iGoldberg Sz Spergel 



t99^ ): IVerde fc Soergell (2002); Cooray & Hu (2000)). In this case there is 



no underlying primordial bispectrum F, but the CMB three-point function 
mi/2"!,2/3"i3 seems to satisfy the factoriz ability condition for a w ide range 
of secondaries (for a survey of shapes see iSmith Sz Zaldarriagal ( 20061 )). 



Summarizing, there is an optimal estimator £ which can be applied to esti- 
mate the amplitude of a three-point signal of specified shape. For the estimator 
to be computationally feasible, the three-point function must satisfy the techni- 
cal requirement of factorizability (Eq. (jlSp ). The factorizability requirement is 
satisfied for the three primordial shapes we will consider in this article (and for 
many secondary shapes as well). Evaluating the estimator also requires a fast 
algorithm for the experiment-specific operation agm {C~^aY"^ which we now 
discuss in more detail in the context of WMAP. 

2.3. Primordial non-Gaussianity in WMAP 

In the previous subsection, the computational problem of evaluating the optimal 
estimator /jy^ was reduced to the "C~^ problem": finding an efficient algorithm 
to compute the inverse (signal-|-noise) weighted map (C~^a)^'" given a harmonic- 
space map aim- 

The operator is experiment-specific. One simplifying feature of WMAP 
is that the instrumental noise can be treated as uncorrelated between pixels 
(more precisely, pixel correlations due to 1/f noise in the detectors are only 
important on large angular scales, where the instrumental noise is much smaller 
than the CMB temperature fluctuations). Even with this simplification, the 



*More precisely, each factorizable term in the primordial bispectrum evolves to a CMB three- 
point function which has an integral representation with factorizable integrad. When the 
integral is approximated by a finite sum, the CMB three-point fun ction will be factorizable in 
the sense defined by Eq. (|15|l . There is an optimization algorithm (jSmith fc Zaldarriaga|[200^ ) 
which can be used to minimize the number of quadrature points needed to approximate the 
integral to a specified level of accuracy. 
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problem is still very challenging for WMAP, due to the large number of pixels 
(A'pix ~ lO'^) needed to represent the WMAP data. A brute-force linear algebra 
approach, representing by a dense iVpix-by-A'^pix matrix, is computationally 
infeasible. Consequently, algorithms for solving the C^^ problem for a large- 
A^pix experiment such as WMAP are based on iterative methods such as the 
conjugate gradient algorithm. 

Iterative algorithms can be used if there is a fast algorithm for applying 
the "forward" operation a — ?> Ca, and a fast algorithm for computing a linear 
operation a — >■ Pa (the "preconditioner" ) such that P approximates as 
closely as possible. The idea behind the iterative algorithms is to compute C~^a 
by iteratively improving a guess x ~ C~^a. In each iteration, the preconditioner 
is applied to the residual vector (a — Cx) to determine a direction in which 
to search for an improved guess. The algorithm terminates when the residual 
vector is sufficiently close to zero. The convergence rate is determined by the 
accuracy of the approximation w P; thus the iterative approach will be 
successful (i.e., fast) if a good preconditioner can be constructed. 

The C^^ problem is an ingredient in many flavors of optimal estimators, 
such as optimal power spectrum estimators, the optimal estimators for primor- 
dial non-Gaussianity in this section, and optimal estimators for gravitational 
lensing to be considered shortly. For this reason, the problem of finding a 
good preconditioner for WMAP has been c onsidered by several au thors (e.g. 
lOh. Spergel &: Hinshawlll999l : lEriksen et al.l[2004: Smith et al.ll2007l 'l. The best 
preconditioner to date was constructed in ISmith et al.l (|2007l ^. using a multigrid 



approach: the preconditioner for the operation is defined by evaluating a 
(partially converged) conjugate gradient operation on a lower-resolution 
version of the WMAP dataset. The preconditioner for the lower-resolution C^^ 
operation is defined using even lower resolution, and so on recursively. Using the 
multigrid preconditioner, the operation a — > C~^a for WMAP can be computed 
in ~15 CPU- minutes. This is sufficiently fast that the optimal estimator f^^ 
can be evaluated in a Monte Carlo pipeline. 

Optima l f'Mr, constraints from WMAP were first reported for the public 5- 
year release ( Smith. Senatore Sz Zaldarriaga|[200^ : ISenatore. Smith &: Zaldarriaga 
I2OICII ). This implementation of the optimal estin iator was subsequent ly incor- 
porated into the WMAP 7- year analysis pipeline (jKomatsu et aLlboid ). At the 
time of this writing, current WMAP constraints on /tvl parameters are: 



f^f' = 32 ±21 
f^f = 26 ±140 
^orthog ^ -202 ±104 (Id error) (17) 



Thus the WMAP data are consistent with Gaussian initial conditions. These 
constraints will improve by a factor ~ 5 in a few years with results from the 
Planck satellite mission. 
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3. Gravitational lensing 
3.1. Introduction 

One of the largest sources of "secondary" CMB anisotropy, or additional anisotropy 
generated after recombination, is gravitational lensing. For purposes of this arti- 
cle, the effect of gravitatio nal lensing can be succinc t ly described as follows. (For 
more detailed reviews, see [Lewis Challinoij(|2006l ):lH anson. Challinor fc LewisI 



( 20091 ).) Gravitational potentials in the late universe deflect photons, so that an 



observer who looks in line-of-sight direction n sees the part of the surface of last 
scattering which lies in direction (n + V0(n)), where the lens potential (/>(n) is 
a 2D field which can be written as a line-of-sight integral. Gravitational lensing 
preserves surface brightness, so that the effect of gravitational lensing is simply 
to move CMB aniostropy around, 

(AT)(n) = (Ar)(n + V(A(n)) (18) 

In this equation and throughout this article, we denote the lensed CMB tem- 
perature by T(n) or a^^, and the unlensed temperature by T(n) or aim- 

Via gravitational lensing, the CMB is indirectly sensitive to the power spec- 
trum Cf'^ of the lens potential. This power spectrum can be written as a line-of- 
sight integral which contains geometric distance factors and the power spectrum 
PiS/{k) of the large-scale Newtonian potential: 



J \ x*x 

In this way, the CMB becomes sensitive to the expansion history and growth 
of structure in the late universe, which adds qualitatively new information. In 
contrast, the unlensed CMB is very sensitive to a single "late universe" paramter, 
the angular diameter distance Z)a(2*) to the redshift of recombination, but is 
otherwise insensitive to the late universe: this is the so-c alled angular diameter 
distance degeneracy ( Zaldarriaga. Spergel &: Seliaklll997l ). CMB lensing breaks 



the angular diameter distance degeneracy and can ultimately provide interesting 
constrai nts on new parameters such as neutrino mass or the dark e nergy equation 
of state (|Stompor Efstathioulll999l : ISmith k Zaldarriagall2006l ). 



How can we reconstruct the lensing power spectrum Cf from observations 
of the CMB? One possibility is to consider the effect of gravitational lensing 
on the CMB temperature power spectrum. Lensing alters Cf^ by smooth- 
ing the acoustic peaks and a dding power in the high-£ tail (ISeliakI Il996l : [13 



20001 : IChallinor Lewisll2005l ) The peak-smoothing effect can be understood as 



degree-scale lenses transferring CMB power between different values of with 
characteristic scale 100. Extra power in the damping tail can be under- 

stood as small-scale lenses generating new anisotropy from a smoothly varying 
CMB backlight. (On small angular scales, the unlensed CMB power spectrum 
is exponentially suppressed but the power spectrum of the lens potential is not). 

A more qualitative effect of CMB lensing is that it generates non-Gaussianity. 
In the next subsection, we will calculate some higher-point correlation functions 
due to lensing. 
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3.2. Three-point and four-point signals from lensing 

Let us first compute the three-point function (TTcj)) between the lens potential 
and two powers of the (lensed) CMB temperature. If we Taylor expand the 
lensed CMB in powers of the lens potential, we get 

a^imi = a<?imi + fhi2t(H2m2^*im\ ^ ^ m ) ^"^^^ 

where the mode-coupling kernel /f^^j^ defined by 

-£l(£l + l)+^2(^2 + l)+-f(^+l) 



V 4^^ V ) ^^^^ 



A short calculation now gives the {TTip) three-point function: 

(a^imia^ama'/'^m) = (/^i^2^C'^^ + /^2<?i<?C'li"^) ^f'^ ^ ^ (22) 

If we replace the lens potential (pim with a different large-scale structure tracer 
field gem, then the correlation between (p and g will generate a {TTg) correlation 
given by: 

a^imiO^amaaL) = (/<?i^2^Q^^ + /^2^i^C'li'^) Cf^ ^ ^ (23) 

We will consider a case where the tracer field is given by number counts of radio 
galaxies in the NVSS survey, in ^3.4.1 

We can also consider the case where the tracer field is the CMB temperature 
itself, where the correlation Cj'^ arises from the ISW effect in a cosmology with 
Qrri 7^ 1- This gives rise to the "ISW-lensing" three-point function, which is 
internal to the CMB: 

{ae,miae^m2ae3ms) = f ht2hCl^ Cj.^ ( ^\ ) + (5 perm.) (24) 

The ISW-lensing signal is the main contribution to the three-point function in- 
duced by gravitational lensing. There is a larger signal in the four-point function. 

mi 0'£2m2 O'isms ^£4m4 ) conn 

(25) 

Ef f r^TT r^TT r^<l>(t> 



X 

m 



V(-l)- ^1 ^3 . .2 t4 . +(11 perm.) 

' \ mi 7713 m I \ m2 m4 —m I ^ ^ ' 



For comparison, forecasts for the Planck satellite show that the ISW-lensing 
three-point signal should be detectable at the ~ 5cj level, whereas the four-point 



10 



Kendrick M. Smith 



signal should be detectable at ~ GOo", and can in fact be used to make a precision 
measurement of Cf"^. 

The three-point and four-point calculations in this section have been some- 
what formal, in the sense that we have written down expressions for the correla- 
tion functions (Eqs. (j22p -(j25 p ) without much interpretation. In the next section 
we will give a more intuitive interpretation for these signals, and also construct 
estimators for extracting them from data. 

3.3. Estimators for gravitational lensing 

The large higher-point signals described in the preceding subsection contain 
information, via the lens potential (j), about the expansion history and growth 
of structure in the late universe. How can we best extract this information, 
i.e. what higher-point estimators should we apply to data to measure these 
signals? We will answer this question for the (TTg) three-point signal (Eq. ([231)), 
the ISW-lensing three-point signal (Eq. (j24p ). and the four-point CMB signal 
(Eq. m)- 

In order to make the analogy with the estimators for primordial non-Gaussianity 
the clearest, let us first consider the ISW-lensing signal. Since this signal is just 
a different shape for the CMB three-point function, the minimum-variance esti- 
mator £ is given by the general form in Eq. (jl3p given previously in the context 
of primordial non-Gaussianity. If we write out the definitions of T and dT for 
the special case of the ISW-lensing three-point function, the estimator becomes: 

c _ ^ \^ f ^TTf~,T<p ( h ^2 h 
^ - p 2^ M^2^3^^2 '-^^ams i mi ma ms 

(.inrii 

^C-^afimr ^(j-l^Y^m2 (C7-la)^3m3 

_^^-1^^2m2,£3m3((^-l^)£imi (26) 

The last two terms in the estimator turn out to be small compared to the first 
two terms, so we will neglect them. Under this approximation, we can rewrite 
£^ in a mathematically equivalent way by introducing the quadratic estimator 
(AA-10), definecdby: 

{N~''(i))im = (-1)"* fhi2tCj^ ( ^ ^ _^ ) 



^We have used the notation M~^(f) rather than <f) because there is a formal sense in which the 
right-hand side of Eq. (|27p is a reconstruction for the lens potential after multiplying by the 
inverse noise covariance matrix of the reconstruction. More precisely, there is an operator 

■^<m,rm' such that {{U'^$)lm) ^ ^ ,e <i> I' -m' aud {{N'^h)lrr,{M'^h)l'-m') =^l-m,l'-m'- In 

the simplified case of all sky coverage and isotropic noise, the operator can be inverted and 
one can write down an estimator (jifm such that {4>im) = 4>tm- In the cut sky case the operator 
A/""^ is noninvertible; for this reason we work with the inverse noise weighted reconstruction 
{M~^4')im in Eq. p7|) and throughout this article. 
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We have arrived at this definition by studying A^-point correlation functions, 
but the quadratic reconstruction N~^(t) was originahy propos ed as a minimum 
variance estimator of the lensing p otential from the CMB fiBernardeaul 1997 



Zaldarriaga &: Seliaklll999l : IHuII2001| ) . Because the lens breaks statistical isotropy 



of the CMB, correlations between different Fourier modes of the CMB can be 
used to estimate the lensing potential. Each Fourier mode (j)^^ of the lensing 
potential can be independently estimated, so the lens reconstruction (f)^^ has 
the degrees of freedom of a map, rather than a scalar quantity. Note that 
our definition of the estimator in Eq. (j27p includes subtraction of the "mean 
field" term ((7~-'^)^i™'i'^2m2_ This can be interpreted as simply subtracting off 
the contribution from spurious anisotropy, induced by anisotropy of the noise 
covariance, when estimating the lensing potential. 

The ISW-lensing estimator £ can then be rewritten: 

^=jT. Cjt'{M~'^)em{C-'a)}^ (28) 

tm 

In this form, the estimator has a simple intuitive interpretation. The ISW- 
lensing signal can be interpreted as a correlation between the (inverse noise 
weighted) lens reconstruction [Af^^c/)) and the (inverse signal-|- noise weighted) 
temperature (C~^a). The optimal estimator is simply the cross power spectrum 
of these two fields, evaluated in a single large bandpower whose "shape" in 
£ is given by the cross spectrum Cj'^. When written in this way (Eq. (j28p ). 
the estimator looks like a two-point quantity, but it is actually a three-point 
statistic in the CMB temperature map, because the lens reconstruction {M~^(j)) 
is quadratic in the CMB. (Note that there is also a one-point counterterm in the 
estimator, coming from the mean field term in the definition of M~^cl).) 

The other higher-point signals from the previous subsection can be treated 
similarly to the ISW-lensing case. Consider next the (TTg) three-point signal 
(Eq. ()23p ). In this case, the optimal estimator is given by: 

^ = I E Cf^{M''^)ern{C-'g)}^ (29) 

Em 

For this estimator to make sense, we must have noisy observations of the large- 
scale structure field gim, with associated noise and signal covariance so that 
{C~^g)im can be defined. Note that the above equation, we have written down 
the optimal estimator for an overall multiple of the fiducial signal. This estimator 
is just the cross power spectrum of the lens reconstruction (pi^ and the observed 
galaxy field, evaluated in a large bandpower with £ weighting given by Cf^. This 
is the appropriate estimator for making a statistical detection of a weak signal, 
but in a case where the (TTg) correlation can be detected with high significance, 
it may be more appropriate to split the estimator in Eq. (|29p into i bands and 

estimate the cross power spectrum Cf^ in independent bandpowers0 



^Let us mention one more subtlety in the TTg estimator (Eq. (|29|l l: inverse noise weighting 
the 4' field is only optimal if the lens reconstruction is noise-dominated {Cf'^ <^ Nf'^). If 
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Finally, we consider the four-point CMB signal given by Eq. (|25p. A nearly- 
optimal estimator is the auto power spectrum of (p, summed over i in one band- 
power with "shape" proportional to the signal power spectrum Cf'^: 

£ = ^J2 cf^ [{Af-'4^)U^f-'^)}m - ((AA-i0),^(AA-i<A)L)] (30) 

£m 

where (•) denotes an expectation value taken over random realizations of the 
CMB. As with the TTg estimator, we can alternately split the i sum into sub- 
ranges and estimate the power spectrum Cf'^ in independent bandpowers. Note 
that this estimator £ is written so that it looks like a two-point quantity, but 
is actually a four-point estimator in the CMB. This estimator is "nearly op- 
timal" because the optimal estimator for Cf"^ contains additional two-point 
counterterms which reduce the variance in the same way that including the 
one-point counterterm in the /nl estimator (Eq. (fTH|) ) reduces the variance, rel- 
ative to an estimator defined w ith only a three-point term. For more details see 



([Smith. Senatore Sz Zaldarriag a 2010). 

In summary, one can construct higher-point estimators for CMB lensing in 
two mathematically equivalent ways. Formally, one can write down expressions 
for three-point or four-point correlation functions generated by lensing, and 
obtain optimal estimators as special cases of general expressions for the optimal 
estimators. More intuitively, one can apply the quadratic lens reconstruction 
estimator tp (Eq. ([27|) ) and then compute cross and auto power spectra. In lens 
reconstruction language, the three-point and four-point correlations generated 
by lensing can be interpreted as either the two-point correlation between cj) and 
another field, or the two-point function of (p itself. 

3.4. Gravitational lensing in WMAP 

In this subsection we will prese n t lens ing results from 3-year WMAP data, orig- 
inally reported in ISmith et al. ( 20071 ). At WMAP resolution, the lens recon- 



struction (f) is highly noise-dominated. Forecasting shows that there is insuffi- 
cient signal-to-noise to detect lensing using the auto power spectrum Cf'^ (i.e. 
the estimator defined in Eq. (j30p l. However, the signal-to-noise can be boosted 
by cross correlating (p with another field which is highly correlated to (p 
less noisy. The best candidate for suc h a cross-correlation is the radio galaxy 
number density from the NVSS survey ( Condon et al.|[l998l ). The NVSS survey 



has large sky coverage (/sky = 0.8), low Poisson noise (A'^gai = 1.8 x 10 ), and 



this assumption is not satisfied then the choice of optimai estimator wiii depend on whether 
optimaiity is defined assuming a fiducial model with Cf"^ — 0. If C^'*|fld is assumed zero (which 
makes sense if one is trying to obtain the most statistically significant detection of a nonzero 
(TTg) signal) then the estimator in Eq. I|29p is optimal. If C^'^|fld is assumed nonzero (which 
makes sense if one is trying to obtain the smallest error bars on Cf^ bandpowers around their 
fiducial values), then the inverse noise weighted field {Af~^4>) should be replaced by a (suitably 

defined) inverse signal+noise weighted field {S + N)~^4> in Eq. ((29}. A similar comment 
applies to the four-point estimator in Eq. (|30|l . This subtlety will be unimportant for the 
WMAP analysis in the next subsection where the lens reconstruction is very noise-dominated. 
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high redshift (2:median = 0.9), making it an excellent match to a large- /sty CMB 
lens reconstruction. 

The main computational problem in implementing the optimal TTg estima- 
tor is the operation appearing in Eq. (j29p . However, a solution to this prob- 
lem has already been described (in the context of primordial non-Gaussianity) 
in El 

Another implementational issue is that the definition of (p given above would 
be computationally prohibitive (the computational cost would be O(^max)) 
evaluated in harmonic space using Eq. (|27|) . However, a fast mathematically 
equivalent expression for (p is given by the following chain of definitions: 

a(n) = ^(C~ia)^^™^y,,„,(n) 

iinii 

A(n) = 5] Q7(C-ia)^^-W,y,,„,(n) 

i^f-'4>)^m = yd'n(V*y4(n))[a(n);3,(n)-(a(n)/3,(n))] (31) 

This is closely analogous to the /nl estimator considered previously: the most 
algebraically straightforward way to write down the estimator (Eq. p3|) ) has 
computational cost ©(^^ax)' ^'^^ there is a faster algorithm based on the specific 
form of the three-point function. 

In Fig. 1, we show the result of evaluating the optimal (TTg) estimator 
(Eq. (j29p ) using CMB measurements from WMAP and galaxy number counts 
from NVSS. We have split the estimator into several £ bands and reported an es- 
timate for Cf^ in each band. The errorbars include systematic errors from num- 
ber density gradients in NVSS, beam effects in WMAP, CMB foregrounds and 
point sources, and thermal SZ. The overall significance (statistical+systematic) 
for detecting nonzero Cf^ is 3Aa. This re sult was the f irst d etection of CMB 
lensing. A detection was also reported bv iHirata et al.l ( 20081 ). using a subop- 



timal estimator but a larger galaxy sample obtained by combining NVSS with 
luminous red galaxies and photometric quasars from SDSS. 



4. Discussion and future prospects 

The theme of this article has been higher-point CMB signals which complement 
the power spectrum. We have studied examples of such signals from inflation- 
ary physics and gravitational lensing in the late universe, constructed optimized 
estimators which are "matched" to various higher-point signals, and reported re- 
sults from WMAP data. Let us conclude with a look to the future, by discussing 
upcoming observational prospects and some of the many unsolved theoretical 
problems in this area. 

For single-field inflation, there is a theorem (Eq. ([5])) which completely char- 
acterizes the most general three-point function which can be generated during 
inflaton. This has not yet been generalized to multifield inflation. The local 
shape f^'i^ is one example of a three-point signal which can arise in multifield 
models but is disallowed in the single-field case, but are there other possibilities? 
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Figure 1. Detection of gravitational lensing in WMAP, obtained by cross- 
correlating the lens reconstruction (j) with galaxy number counts from NVSS 
using the optimal estimator (Eq. (1291) ). The estimated power spectrum Cf^ is 
consistent with the fiducial model in the WMAP3 cosmology. The two sets of 
error bars represent statistical and (statistical + systematic) errors. Sources 
of systematic error include number density gradients in NVSS, beam effects in 
WMAP, Galactic microwave foregrounds, resolved and unresolved CMB point 
sources, and the thermal Sunyaev-Zeldovich effect. The significance (statisti- 
cal+systematic) of the detection is 3.4(T, by fitting the estimated bandpowers 
to an overall multiple of the fiducial power spectrum Cf^. (From [Smith et aD 
(|2M3).) 
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Is it always sufficient to look for primordial three-point signals, or do there exist 
inflationary models which generate detectably large four-point signals with no 
accompanying detectable three-point signal? 

So far we have not mentioned CMB polarization. Future generations of low- 
noise polarization experiments will be exquisite probes of CMB secondaries such 

as lensing, since first-order perturbative effects with scalar sources only generate 

an E- mode in polarization ( Seliak Sz Zaldarriaga|[l997l : lKamionkowski. Kosowskv &: StebbinsI 
19971 ). whereas secondary effects tend to generate a mixture of E and B-modes. 



For example, in the limit of low instrumental noise, the reconstructed lens po- 
tential 4> from CMB polarization extends t o much smaller angula r scales than 
would be possible using CMB temperature ( Okamoto Sz Hu 20031 ) . As another 
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example, secondary polarization generated by Thomson scattering of CMB pho- 
tons by HIT bubbles during the epoch of inhomogeneous reionization is a mixture 
of E and B-modes and generates a higher-point signal whi ch can be extracted 
by suitably constructed estimators ( Dvorkin Sz Smithll2009l ). 

At the time of this writing, CMB lensing has been detected in cross-correlation 
with large-scale structure, with low statistical significance (3.4(t). The situation 
will change dramatically in a few years with lensing results from Planck and 
ground-based experiments like SPT and ACT, which can constrain lensing "in- 
ternally" (i.e. without a cross-correlation tracer) at the few percent level and 
obtain interesting constraints on parameters such as neutrino mass. The CMB 
lens reconstruction (p from Planck will be the first all-sky lensing map with in- 
teresting signal-to-noise, and will measure a redshift range which is difficult to 
measure with other probes of lensing such as cosmic shear. Gravitational lensing 
will soon be an indirect but important scientific product of experiments which 
measure the small-scale CMB, in much the same way that weak lensing (via 
galaxy ellipticities) is an important product of wide-field optical surveys like 
SDSS. 

On small angular scales, secondary CMB anisotropy generated well after 
recombination dominates the primary CMB. In this regime, the CMB should 
be thought of in a different way: each source of secondary anisotropy is a non- 
Gaussian signal whose statistical distribution is different from the other sec- 
ondaries. Exploring this new observational frontier will require new statistical 
tools; the higher-point estimators presented here represent one approach to this 
problem. Ideally, we would like to have a complete set of estimators which can 
separate the various secondaries (lensing, inhomogeneous reionization, kinetic 
SZ, etc.) from each other, but such a framework has not yet been developed. 
There are many unsolved problems in the field, but the potential scientific re- 
turns from new measurements of CMB temperature and polarization in the next 
few years are very exciting. 
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